clear; close all; clc;

% 设置迭代次数和去除瞬态
N = 5000;
transient = 1000;

%% 情况 (a)：混沌爆发，参数 (alpha, sigma, mu) = (4.5, 0.14, 0.001)
alpha1 = 4.5; 
sigma1 = 0.14; 
mu1 = 0.001;

x1 = zeros(1, N); 
y1 = zeros(1, N);

% 初始值
x1(1) = 0.01; 
y1(1) = 0;

for n = 1 : N-1
    % Rulkov映射中 F(x, y) 的分段定义
    if x1(n) <= 0
        F = alpha1/(1 - x1(n)) + y1(n);
    elseif x1(n) < alpha1 + y1(n)
        F = alpha1 + y1(n);
    else
        F = -1;
    end
    
    x1(n+1) = F;
    y1(n+1) = y1(n) - mu1*(x1(n) + 1) + mu1*sigma1;
end

%% 情况 (b)：周期爆发，参数 (alpha, sigma, mu) = (6, -0.1, 0.001)
alpha2 = 6; 
sigma2 = -0.1; 
mu2 = 0.001;

x2 = zeros(1, N); 
y2 = zeros(1, N);

% 初始值
x2(1) = 0.01; 
y2(1) = 0;

for n = 1 : N-1
    if x2(n) <= 0
        F = alpha2/(1 - x2(n)) + y2(n);
    elseif x2(n) < alpha2 + y2(n)
        F = alpha2 + y2(n);
    else
        F = -1;
    end
    
    x2(n+1) = F;
    y2(n+1) = y2(n) - mu2*(x2(n) + 1) + mu2*sigma2;
end

%% 绘图：忽略初始瞬态，绘制两种放电模式
figure;

% (a) 绘制区间 [2500, 3000] 混沌爆发行为
subplot(2,1,1);
plot(x1, 'b');            % 先完整绘制序列
xlim([2500 3500]);        % 限制横坐标区间
xlabel('迭代次数 n ');
ylabel('x');
title('(α, σ, μ) = (4.5, 0.14, 0.001)下混沌爆发行为');

% (b) 绘制区间 [2500, 3500] 周期爆发行为
subplot(2,1,2);
plot(x2, 'r');
xlim([2500 3500]);
xlabel('迭代次数 n ');
ylabel('x');
title('(α, σ, μ) = (6, -0.1, 0.001)下周期爆发行为');

